MSK SPECTRUM data freeze: major cell type embeddings
Florian Uhlitz
2021-03-09
library(magrittr)
library(tidyverse)
library(Seurat)
library(readxl)
library(cowplot)
library(colorblindr)
library(viridis)
library(magick, lib.loc = "/home/uhlitzf/miniconda3/lib/R/library")
library(ggpubr)## load global vars:
source("_src/global_vars.R")
# meta_tbl
# clrs
# markers_v7
# markers_v7_super
# cell_type_super_lookup## load full seurat objects with expression data
# seu_obj_tc <- read_rds("/work/shah/uhlitzf/data/SPECTRUM/freeze/v5/T.cell_processed_filtered_sub.rds")
# seu_obj_cc <- read_rds(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v5/Ovarian.cancer.cell_processed_filtered.rds"))
seu_obj_cc <- read_rds("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/Ovarian.cancer.super_processed_filtered_annotated.rds")
marker_tbl <- read_tsv("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/supplementary_tables/Ovarian.cancer.super_marker_table_annotated.tsv")
marker_tbl_top <- marker_tbl %>%
filter(avg_logFC > 0.5,
p_val_adj < 0.01,
pct.1 > 0.2,
pct.2 < 0.8,
!is.na(cluster_label)) %>%
group_by(cluster_label) %>%
slice(1:50)
myfeatures <- c("UMAP_1", "UMAP_2", "umapharmony_1", "umapharmony_2", "sample", "doublet", "nCount_RNA", "nFeature_RNA", "percent.mt", "doublet_score", "cell_type")
my_subtypes <- names(clrs$cluster_label$Ovarian.cancer.super)
coi <- "Ovarian.cancer.super"plot_data <- cbind(cell_id = colnames(seu_obj_cc), FetchData(seu_obj_cc, c("umapharmony_1", "umapharmony_2", "umappca_1", "umappca_2", "RNA_snn_res.0.2", "sample", "cluster_label", grep("pathway", colnames(seu_obj_cc@meta.data), value = T)))) %>%
as_tibble() %>%
left_join(meta_tbl, by = "sample") %>%
filter(!is.na(consensus_signature)) %>%
mutate(cluster_label = ordered(cluster_label, levels = names(clrs$cluster_label$Ovarian.cancer.super)),
cluster_number = as.numeric(cluster_label))
base_umap <- ggplot(plot_data) +
coord_fixed() +
NoAxes() +
theme(legend.position = c(0, 1),
legend.justification = c("left", "top"),
legend.box.just = "left",
legend.margin = margin(1, 1, 1, 1),
#panel.border = element_rect(linetype = 1, color = "black", size = 1),
legend.text = element_text(size = 14, margin = margin(0, 10, 0, 0)),
legend.spacing.x = unit(0, "npc"),
legend.title = element_blank(),
plot.title = element_text(hjust = 0.5, vjust = 0.5, face = "plain", size = 22))
arrow <- arrow(angle = 20, type = "closed", length = unit(0.1, "npc"))
umap_coord_anno <- ggplot(tibble(group = c("UMAP1", "UMAP2"),
x = c(0, 0), xend = c(1, 0),
y = c(0, 0), yend = c(0, 1),
lx = c(0.5, -0.15), ly = c(-0.15, 0.5),
angle = c(0, 90))) +
geom_segment(aes(x, y, xend = xend, yend = yend, group = group),
arrow = arrow, size = 1, lineend = "round") +
geom_text(aes(lx, ly, label = group, angle = angle), size = 4) +
theme_void() +
coord_fixed(xlim = c(-0.3, 1), ylim = c(-0.3, 1))
add_umap_coord <- function(gg_obj) {
p <- ggdraw() +
draw_plot(gg_obj, x = 0, y = 0, width = 1, height = 1) +
draw_plot(umap_coord_anno, x = -0.015, y = -0.02, width = 0.2, height = 0.2)
return(p)
}pt.size <- 0.1
pt.size2 <- 0.2
pt.size.mini <- 0.01
pt.alpha <- 0.05
pt.alpha.mini <- 0.02
median_tbl <- plot_data %>%
group_by(patient_id_short) %>%
summarise(umappca_1 = median(umappca_1),
umappca_2 = median(umappca_2))
median_tbl_cluster <- plot_data %>%
group_by(cluster_label, cluster_number) %>%
summarise(umapharmony_1 = median(umapharmony_1),
umapharmony_2 = median(umapharmony_2))
umap_pca_mutsig <- base_umap +
geom_point(aes(umappca_1, umappca_2, color = consensus_signature),
size = pt.size, alpha = pt.alpha) +
# geom_text(aes(umappca_1, umappca_2, label = patient_id_short), data = median_tbl) +
geom_label(aes(umappca_1, umappca_2, label = patient_id_short), color = "black",
data = median_tbl, label.size = unit(0, "mm"), label.r = unit(0, "mm"),
alpha = 0.5) +
scale_color_manual(values = clrs$consensus_signature) +
guides(color = guide_legend(override.aes = list(size = 4, alpha = 1),
ncol = 1,
label.position = "right")) +
labs(title = "Uncorrected")
umap_mutsig <- base_umap +
geom_point(aes(umapharmony_1, umapharmony_2, color = consensus_signature),
size = pt.size2, alpha = pt.alpha) +
scale_color_manual(values = clrs$consensus_signature) +
guides(color = guide_legend(override.aes = list(size = 4, alpha = 1),
ncol = 1,
label.position = "right")) +
guides(color = F) +
labs(title = "Corrected (harmony)")
umap_pca_cluster <- base_umap +
geom_point(aes(umappca_1, umappca_2, color = cluster_label),
size = pt.size, alpha = pt.alpha) +
geom_text(aes(umappca_1, umappca_2, label = patient_id_short), data = median_tbl) +
scale_color_manual(values = clrs$cluster_label$Ovarian.cancer.super) +
guides(color = guide_legend(override.aes = list(size = 4, alpha = 1),
ncol = 1,
label.position = "right")) +
labs(title = "Cluster")
umap_cluster <- base_umap +
geom_point(aes(umapharmony_1, umapharmony_2, color = cluster_label),
size = pt.size2, alpha = pt.alpha) +
scale_color_manual(values = clrs$cluster_label$Ovarian.cancer.super) +
geom_point(aes(umapharmony_1, umapharmony_2),
color = "white", alpha = 0.5, size = 6,
data = median_tbl_cluster) +
geom_text(aes(umapharmony_1, umapharmony_2, label = cluster_number),
color = "black",
data = median_tbl_cluster) +
guides(color = guide_legend(override.aes = list(size = 4, alpha = 1),
ncol = 1,
label.position = "right")) +
guides(color = F) +
labs(title = "Cluster")
cluster_legend <- cowplot::get_legend(umap_pca_cluster)
umap_pca_jak_stat <- base_umap +
# geom_point(aes(umappca_1, umappca_2), color = "grey80",
# size = pt.size, alpha = pt.alpha,
# data = filter(plot_data, JAK.STAT.pathway <= 0)) +
geom_point(aes(umappca_1, umappca_2, color = JAK.STAT.pathway),
size = pt.size, alpha = pt.alpha,
data = mutate(plot_data, JAK.STAT.pathway = ifelse(JAK.STAT.pathway > 4, 4, JAK.STAT.pathway))) +
geom_text(aes(umappca_1, umappca_2, label = patient_id_short), data = median_tbl) +
scale_color_gradientn(colours = viridis(9), breaks = c(0, 2, 4), limits = c(min(plot_data$JAK.STAT.pathway), 4), labels = c(0, 2, "≥4")) +
labs(title = "JAK-STAT signaling")
umap_jak_stat <- base_umap +
# geom_point(aes(umapharmony_1, umapharmony_2), color = "grey80",
# size = pt.size, alpha = pt.alpha,
# data = filter(plot_data, JAK.STAT.pathway <= 0)) +
geom_point(aes(umapharmony_1, umapharmony_2, color = JAK.STAT.pathway),
size = pt.size2, alpha = pt.alpha,
data = mutate(plot_data, JAK.STAT.pathway = ifelse(JAK.STAT.pathway > 4, 4, JAK.STAT.pathway))) +
scale_color_gradientn(colours = viridis(9), breaks = c(0, 2, 4), limits = c(min(plot_data$JAK.STAT.pathway), 4), labels = c(0, 2, "≥4")) +
labs(title = "JAK-STAT signaling")source("_src/comp_plot.R")
cluster_comp <- plot_data %>%
filter(consensus_signature != "Undetermined") %>%
mutate(sort_short_x = str_replace_all(sort_short, "U", "CD45-")) %>%
mutate(sample_id = sample) %>%
group_by(cluster_label, consensus_signature, sample_id, sort_short_x, tumor_supersite, therapy) %>%
tally %>%
group_by(consensus_signature, sample_id, sort_short_x, tumor_supersite, therapy) %>%
mutate(nrel = n/sum(n)*100) %>%
ungroup %>%
mutate(cell_type = "Ovarian.cancer.cell")
plist1 <- default_comp_grid_list(filter(cluster_comp, sort_short_x == "CD45-"),
cluster_label, "Cancer.cell.1", cluster_label,
vec_plot = F, site_box = T,
super_type = "Ovarian.cancer.super")
plist2 <- default_comp_grid_list(filter(cluster_comp, sort_short_x == "CD45-"),
cluster_label, "Cancer.cell.2", cluster_label,
vec_plot = F, site_box = T, yaxis = F,
super_type = "Ovarian.cancer.super")
plist3 <- default_comp_grid_list(filter(cluster_comp, sort_short_x == "CD45-"),
cluster_label, "Cancer.cell.3", cluster_label,
vec_plot = F, site_box = T, yaxis = F,
super_type = "Ovarian.cancer.super")
plist4 <- default_comp_grid_list(filter(cluster_comp, sort_short_x == "CD45-"),
cluster_label, "Cancer.cell.4", cluster_label,
vec_plot = F, site_box = T, yaxis = T,
super_type = "Ovarian.cancer.super")
plist5 <- default_comp_grid_list(filter(cluster_comp, sort_short_x == "CD45-"),
cluster_label, "Cancer.cell.5", cluster_label,
vec_plot = F, site_box = T, yaxis = F,
super_type = "Ovarian.cancer.super")
plist6 <- default_comp_grid_list(filter(cluster_comp, sort_short_x == "CD45-"),
cluster_label, "Cancer.cell.6", cluster_label,
vec_plot = F, site_box = T, yaxis = F,
super_type = "Ovarian.cancer.super")
pcomp_grid_p1 <- plot_grid(plotlist = plist1,
ncol = 1, align = "v",
rel_heights = c(0.2, 0.2, 0.22, 0.38))
pcomp_grid_p2 <- plot_grid(plotlist = plist2,
ncol = 1, align = "v",
rel_heights = c(0.2, 0.2, 0.22, 0.38))
pcomp_grid_p3 <- plot_grid(plotlist = plist3,
ncol = 1, align = "v",
rel_heights = c(0.2, 0.2, 0.22, 0.38))
pcomp_grid_p4 <- plot_grid(plotlist = plist4,
ncol = 1, align = "v",
rel_heights = c(0.2, 0.2, 0.22, 0.38))
pcomp_grid_p5 <- plot_grid(plotlist = plist5,
ncol = 1, align = "v",
rel_heights = c(0.2, 0.2, 0.22, 0.38))
pcomp_grid_p6 <- plot_grid(plotlist = plist6,
ncol = 1, align = "v",
rel_heights = c(0.2, 0.2, 0.22, 0.38))
pcomp_grid_full_1 <- plot_grid(pcomp_grid_p1, pcomp_grid_p2, ggdraw(),
pcomp_grid_p3, ggdraw(),
nrow = 1,
rel_widths = c(0.48, 0.25, 0.01, 0.25, 0.01))
pcomp_grid_full_2 <- plot_grid(pcomp_grid_p4, pcomp_grid_p5, ggdraw(),
pcomp_grid_p6, ggdraw(),
nrow = 1,
rel_widths = c(0.48, 0.25, 0.01, 0.25, 0.01))
pcomp_grid_full_1
pcomp_grid_full_2
ggsave("_fig/003_cancer_cell/003_comp_plot_1.pdf", pcomp_grid_full_1, width = 6, height = 5.5)
ggsave("_fig/003_cancer_cell/003_comp_plot_2.pdf", pcomp_grid_full_2, width = 6, height = 5.5)plist3_2 <- default_comp_grid_list(filter(cluster_comp, sort_short_x == "CD45-"),
cluster_label, "Cancer.cell.3", cluster_label,
vec_plot = F, site_box = T, yaxis = T,
super_type = "Ovarian.cancer.super")
pcomp_grid_p3_2 <- plot_grid(plotlist = plist3_2,
ncol = 1, align = "v",
rel_heights = c(0.2, 0.2, 0.22, 0.38))
pcomp_grid_p3_2
ggsave("_fig/003_cancer_cell/003_comp_plot_cluster3.pdf", width = 3, height = 5.5)cancer_grid_pca <- ggdraw() +
draw_plot(add_umap_coord(umap_pca_mutsig),
x = 0, y = 0, width = 0.25, height = 1) +
draw_plot(add_umap_coord(umap_pca_cluster + guides(color = F)),
x = 0.24, y = 0, width = 0.25, height = 1) +
draw_grob(cluster_legend, x = 0.5, y = -0.15, height = 1) +
draw_plot(add_umap_coord(umap_pca_jak_stat),
x = 0.7, y = 0, width = 0.25, height = 1)
cancer_grid_harmony <- ggdraw() +
draw_plot(add_umap_coord(umap_mutsig),
x = 0, y = 0, width = 0.25, height = 1) +
draw_plot(add_umap_coord(umap_cluster + guides(color = F)),
x = 0.24, y = 0, width = 0.25, height = 1) +
draw_plot(pcomp_grid_p3_2,
x = 0.5, y = 0, width = 0.17, height = 1.02) +
draw_plot(add_umap_coord(umap_jak_stat),
x = 0.7, y = 0, width = 0.25, height = 1)
cancer_grid <- ggdraw() +
draw_plot(add_umap_coord(umap_pca_mutsig),
x = 0, y = 0, width = 0.25, height = 1) +
draw_plot(add_umap_coord(umap_mutsig + guides(color = F)),
x = 0.24, y = 0, width = 0.25, height = 1) +
draw_plot(add_umap_coord(umap_cluster + guides(color = F)),
x = 0.48, y = 0, width = 0.25, height = 1) +
draw_plot(add_umap_coord(umap_jak_stat),
x = 0.73, y = 0, width = 0.25, height = 1)
cancer_grid_pca
cancer_grid_harmony
cancer_grid
# ggsave("_fig/003_cancer_cell/003_umap_grid.pdf", cancer_grid, width = 20, height = 5)
ggsave("_fig/003_cancer_cell/003_umap_grid_pca.png", cancer_grid_pca, width = 20, height = 5)
ggsave("_fig/003_cancer_cell/003_umap_grid_harmony.png", cancer_grid_harmony, width = 20, height = 5)
ggsave("_fig/003_cancer_cell/003_umap_grid.png", cancer_grid, width = 20, height = 5)1 Cluster marker heatmap
plot_data_markers <- as_tibble(FetchData(seu_obj_cc, c("cluster_label", myfeatures, unique(marker_tbl_top$gene)))) %>%
gather(gene, value, -c(1:(length(myfeatures)+1))) %>%
left_join(meta_tbl, by = "sample") %>%
mutate(cluster_label = ordered(cluster_label, levels = names(clrs$cluster_label$Ovarian.cancer.super))) %>%
group_by(cluster_label, gene) %>%
summarise(value = mean(value, na.rm = T)) %>%
group_by(gene) %>%
mutate(value = scales::rescale(value)) %>%
left_join(select(marker_tbl_top, cluster_label_x = cluster_label, gene), by = "gene") %>%
## reverse names vector to flip row-order in heatmap
mutate(cluster_label_x = ordered(cluster_label_x, levels = names(clrs$cluster_label$Ovarian.cancer.super))) %>%
na.omit()library(ComplexHeatmap)
highlight_genes <- marker_tbl_top %>%
group_by(cluster_label) %>%
slice(1:2) %>%
## reverse levels also here for row-order flip
mutate(cluster_label_x = ordered(cluster_label, levels = names(clrs$cluster_label[[coi]]))) %>%
ungroup() %>%
select(cluster_label_x, gene) %>%
na.omit %>%
mutate(highlight = T)
plot_data_markers_mat <- plot_data_markers %>%
spread(cluster_label, value) %>%
left_join(highlight_genes, by = c("gene", "cluster_label_x")) %>%
ungroup %>%
arrange(desc(cluster_label_x), gene) %>%
select(cluster_label_x, gene, highlight, everything())
ha_row <- rowAnnotation(
`Cell type` = plot_data_markers_mat$cluster_label_x,
col = list(`Cell type` = clrs$cluster_label[[coi]]),
show_legend = F,
annotation_name_side = "top"
)
cluster_names <- colnames(plot_data_markers_mat)[-c(1:3)]
ha_col <- columnAnnotation(
`Cell type` = anno_block(gp = gpar(fill = clrs$cluster_label[[coi]], col = NA),
labels = as.character(1:10)),
show_legend = F,
annotation_name_side = "left",
annotation_name_rot = 0
)
gene_idx <- plot_data_markers_mat$highlight == TRUE
ha_genes <- rowAnnotation(
link = anno_mark(
at = which(gene_idx),
labels = plot_data_markers_mat$gene[which(gene_idx)],
labels_gp = gpar(fontsize = 10), padding = unit(1, "mm"),
side = "left",
labels_rot = 0
)
)
marker_heatmap <- Heatmap(
as.matrix(plot_data_markers_mat[,-c(1:3)]),
heatmap_legend_param = list(
title = "Scaled expression",
title_position = "leftcenter-rot"
),
row_order = 1:length(plot_data_markers_mat$cluster_label_x),
row_split = plot_data_markers_mat$cluster_label_x,
column_split = 1:length(colnames(plot_data_markers_mat)[-c(1:3)]),
column_order = colnames(plot_data_markers_mat)[-c(1:3)],
column_names_side = "top",
right_annotation = ha_row,
left_annotation = ha_genes,
top_annotation = ha_col,
cluster_rows = F,
row_title = NULL,
column_title = NULL,
col = viridis(9)
)
# marker_heatmap
heatmap_grob <- grid.grabExpr(draw(marker_heatmap), width = 3.75, height = 5)
heat_grid <- ggdraw() +
draw_grob(heatmap_grob)
heat_grid
ggsave("_fig/003_cancer_cell/003_cluster_marker_heatmap.pdf", heat_grid, width = 3.75, height = 5)2 Differential pathway expression in cancer cells
Are there certain pathways activated in cancer cells…
- of patients with high vs low dysfunctional T cells?
- of HRD vs FBI patients
set.seed(42)
sampled_cell_ids <- sample(colnames(seu_obj_cc), 10000)
seu_obj_cc_sub <- subset(seu_obj_cc, cells = sampled_cell_ids)
plot_data <- FetchData(seu_obj_cc, c("umapharmony_1", "umapharmony_2", "sample", "cluster_label", grep("pathway|module", colnames(seu_obj_cc@meta.data), value = T))) %>%
as_tibble() %>%
gather(pathway, score, -c(1:4)) %>%
left_join(meta_tbl, by = "sample") %>%
filter(sort_short == "CD45-", therapy == "pre-Rx") %>%
mutate(pathway = str_remove_all(pathway, "\\.pathway")) %>%
mutate(cluster_label = str_remove_all(cluster_label, "\\.cell")) %>%
filter(consensus_signature != "Undetermined")
cut_value <- 2
pathway_summary_wrapper <- . %>%
summarise(mean_score = mean(score),
median_score = median(score)) %>%
mutate(median_cut = ifelse(median_score > cut_value, cut_value,
ifelse(median_score < -cut_value, -cut_value,
median_score))) %>%
mutate(mean_cut = ifelse(mean_score > cut_value, cut_value,
ifelse(mean_score < -cut_value, -cut_value,
mean_score)))
# plot_data_summary_patient <- plot_data %>%
# group_by(patient_id_short, consensus_signature, cluster_label,
# pathway, tumor_supersite) %>%
# pathway_summary_wrapper
#
# plot_data_summary_mutsig <- plot_data %>%
# group_by(consensus_signature, pathway) %>%
# pathway_summary_wrapper
#
# plot_data_summary_mutsig_cluster <- plot_data %>%
# group_by(consensus_signature, cluster_label, pathway) %>%
# pathway_summary_wrapper
plot_data_summary_cluster <- plot_data %>%
group_by(cluster_label, pathway) %>%
pathway_summary_wrapper
plot_data_summary_mutsig_patient <- plot_data %>%
group_by(consensus_signature, patient_id_short, pathway) %>%
pathway_summary_wrapper
plot_data_summary_mutsig_sample <- plot_data %>%
group_by(consensus_signature, sample, pathway) %>%
pathway_summary_wrappercommon_heat_layers <- list(
scale_fill_gradient2(low = scales::muted("blue"), high = scales::muted("red"),
na.value = "grey10",
breaks = c(-cut_value, 0, cut_value),
labels = c(paste0("≤-", cut_value), "0", paste0("≥", cut_value)),
limits = c(-cut_value, cut_value)),
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5),
axis.ticks = element_blank(),
axis.line = element_blank(),
axis.title = element_blank(),
plot.margin = margin(0, 0, 0, 0))
)
# ggplot(plot_data_summary_patient) +
# geom_tile(aes(cluster_label, patient_id_short, fill = median_cut)) +
# facet_grid(consensus_signature~pathway, scales = "free", space = "free") +
# common_heat_layers
#
# ggplot(plot_data_summary_mutsig_cluster) +
# geom_tile(aes(pathway, cluster_label, fill = median_cut)) +
# facet_grid(~consensus_signature, scales = "free", space = "free") +
# common_heat_layers
#
# ggplot(plot_data_summary_mutsig) +
# geom_tile(aes(consensus_signature, pathway, fill = median_score)) +
# common_heat_layers
pw_plot1 <- plot_data_summary_cluster %>%
mutate(facet_helper = "") %>%
filter(!str_detect(pathway, "module")) %>%
ggplot() +
geom_tile(aes(cluster_label, pathway, fill = mean_score)) +
common_heat_layers +
facet_grid(~facet_helper, scales = "free", space = "free") +
labs(x = "Cluster", y = "Pathway", fill = "Mean\nPROGENy\nscore") +
scale_fill_gradient2(low = scales::muted("blue"), high = scales::muted("red"),
# labels = c(paste0("≤-", cut_value), "0", paste0("≥", cut_value)),
# breaks = c(-cut_value, 0, cut_value),
na.value = "grey10") +
theme(axis.text.x = element_blank())
pw_plot1_anno <- plot_data_summary_cluster %>%
mutate(facet_helper = "") %>%
filter(!str_detect(pathway, "module")) %>%
ggplot() +
geom_tile(aes(cluster_label, facet_helper, fill = cluster_label)) +
common_heat_layers +
scale_fill_manual(values = clrs$cluster_label$Ovarian.cancer.super %>% setNames(str_remove_all(names(.), "\\.cell"))) +
facet_grid(~facet_helper, scales = "free", space = "free") +
theme(axis.text.y = element_blank(),
strip.text = element_blank()) +
guides(fill = F)
# pw_plot2 <- ggplot(plot_data_summary_mutsig_patient) +
# geom_tile(aes(patient_id_short, pathway, fill = mean_cut)) +
# facet_grid(~consensus_signature, scales = "free", space = "free") +
# common_heat_layers +
# labs(x = "Patient", y = "", fill = "Mean\nPROGENy\nscore") +
# theme(axis.text.x = element_blank())
# pw_plot2_anno <- ggplot(mutate(plot_data_summary_mutsig_patient, facet_helper = "")) +
# geom_tile(aes(patient_id_short, facet_helper, fill = consensus_signature)) +
# common_heat_layers +
# scale_fill_manual(values = clrs$consensus_signature) +
# facet_grid(~consensus_signature, scales = "free", space = "free") +
# theme(axis.text.y = element_blank(),
# strip.text = element_blank()) +
# guides(fill = F)
pw_grid_left <- plot_grid(pw_plot1, pw_plot1_anno, ggdraw(),
ncol = 1, align = "v", axis = "lrtb",
rel_heights = c(0.65, 0.35, 0))
# pw_grid_right <- plot_grid(pw_plot2, pw_plot2_anno, ggdraw(),
# ncol = 1, align = "v", axis = "lrtb",
# rel_heights = c(0.7, 0.1, 0.2))
comparison_data_jakstat <- filter(plot_data_summary_mutsig_sample, pathway == "JAK.STAT") %>%
rename(score = mean_score) %>%
ungroup
# comparison_data_jakstat <- filter(plot_data_summary_mutsig_sample, pathway == "JAK.STAT") %>% rename(score = mean_score) %>% ungroup %>% mutate(consensus_signature = as.numeric(consensus_signature))
#
# test_result_jakstat <- compare_means(
# consensus_signature ~ score, ref.group = ".all.",
# p.adjust.method = "fdr", method = "wilcox.test",
# data = comparison_data_jakstat
# ) %>%
# mutate(y.position = 3.5)
pw_boxplot_mutsig_jakstat <- filter(plot_data, pathway == "JAK.STAT") %>%
ggplot(aes(consensus_signature, score)) +
geom_violin(aes(consensus_signature, score, fill = consensus_signature), color = "white", adjust = 2, alpha = 0.5, width = 1.5) +
geom_boxplot(aes(consensus_signature, score, color = consensus_signature),
width = 0.5, size = 0.75, outlier.shape = NA) +
geom_boxplot(aes(consensus_signature, score, fill = consensus_signature),
color = "white", width = 0.5, outlier.shape = NA, size = 0.5) +
stat_compare_means(aes(consensus_signature, score, color = "red"),
ref.group = ".all.", data = comparison_data_jakstat,
label = "p.signif", label.y = 3.5, hide.ns = T) +
# stat_pvalue_manual(test_result_jakstat, label = "p.adj") +
stat_compare_means(aes(consensus_signature, score),
data = comparison_data_jakstat,
label.y = 4.3, label.x = 1.5, label.sep = "\n") +
facet_wrap(~pathway) +
scale_color_manual(values = clrs$consensus_signature) +
scale_fill_manual(values = clrs$consensus_signature) +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5)) +
labs(x = "", y = "PROGENy score") +
coord_cartesian(ylim = c(-1, 5))
comparison_data_module <- filter(plot_data_summary_mutsig_sample, pathway %in% c("IFNg.signaling.module", "ISG.module")) %>%
mutate(pathway = c(`IFNg.signaling.module` = "IFNg", `ISG.module` = "IFNa")[pathway]) %>%
rename(score = mean_score)
pw_boxplot_mutsig_module <- filter(plot_data, pathway %in% c("IFNg.signaling.module", "ISG.module")) %>%
mutate(pathway = c(`IFNg.signaling.module` = "IFNg", `ISG.module` = "IFNa")[pathway]) %>%
ggplot(aes(consensus_signature, score)) +
geom_violin(aes(consensus_signature, score, fill = consensus_signature), color = "white", adjust = 2, alpha = 0.5, width = 1.5) +
geom_boxplot(aes(consensus_signature, score, color = consensus_signature),
width = 0.5, size = 0.75, outlier.shape = NA) +
geom_boxplot(aes(consensus_signature, score, fill = consensus_signature),
color = "white", width = 0.5, outlier.shape = NA, size = 0.5) +
stat_compare_means(aes(consensus_signature, score),
ref.group = ".all.", data = comparison_data_module,
label = "p.signif", label.y = 0.9, hide.ns = T) +
stat_compare_means(aes(consensus_signature, score),
data = comparison_data_module,
label.y = 1.1, label.x = 1.5, label.sep = "\n\n") +
facet_wrap(~pathway) +
scale_color_manual(values = clrs$consensus_signature) +
scale_fill_manual(values = clrs$consensus_signature) +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5)) +
labs(x = "", y = "Module score") +
coord_cartesian(ylim = c(-0.5, 1.2))
plot_box_grid <- plot_grid(pw_boxplot_mutsig_jakstat + remove_guides,
pw_boxplot_mutsig_module + remove_guides,
nrow = 1, align = "h",
rel_widths = c(c(0.37, 0.63)))
pw_grid_full <- ggdraw() +
draw_plot(pw_grid_left, x = 0.01, y = 0, width = 0.41, height = 1) +
draw_plot(plot_box_grid, x = 0.46, y = 0.05, width = 0.5, height = 0.9)
pw_grid_full
ggsave("_fig/003_cancer_cell/003_pathway_heatmap_boxplot.pdf", pw_grid_full,
width = 10, height = 5)plot_data_patient_jakstat <- filter(plot_data, pathway == "JAK.STAT") %>%
mutate(patient_id_short = ordered(patient_id_short, levels = unique(arrange(plot_data_summary_mutsig_patient, pathway != "JAK.STAT", median_score)$patient_id_short)))
pw_boxplot_patient_jakstat <- ggplot(plot_data_patient_jakstat) +
geom_violin(aes(patient_id_short, score, fill = consensus_signature), color = "white", adjust = 2, alpha = 0.5, width = 1.5) +
geom_boxplot(aes(patient_id_short, score, color = consensus_signature),
width = 0.5, size = 0.75, outlier.shape = NA) +
geom_boxplot(aes(patient_id_short, score, fill = consensus_signature),
color = "white", width = 0.5, outlier.shape = NA, size = 0.5) +
facet_wrap(~pathway) +
scale_color_manual(values = clrs$consensus_signature) +
scale_fill_manual(values = clrs$consensus_signature) +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5)) +
labs(x = "Patient", y = "PROGENy score",
fill = "Mutational\nsignature",
color = "Mutational\nsignature") +
ylim(c(-1, 5))
pw_boxplot_patient_module <- filter(plot_data, pathway %in% c("IFNg.signaling.module", "ISG.module")) %>%
mutate(pathway = c(`IFNg.signaling.module` = "IFNg", `ISG.module` = "IFNa")[pathway],
patient_id_short = ordered(patient_id_short, levels = levels(plot_data_patient_jakstat$patient_id_short))) %>%
ggplot() +
geom_violin(aes(patient_id_short, score, fill = consensus_signature), color = "white", adjust = 2, alpha = 0.5, width = 1.5) +
geom_boxplot(aes(patient_id_short, score, color = consensus_signature),
width = 0.5, size = 0.75, outlier.shape = NA) +
geom_boxplot(aes(patient_id_short, score, fill = consensus_signature),
color = "white", width = 0.5, outlier.shape = NA, size = 0.5) +
facet_wrap(~pathway, ncol = 1) +
scale_color_manual(values = clrs$consensus_signature) +
scale_fill_manual(values = clrs$consensus_signature) +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5)) +
labs(x = "Patient", y = "Module score",
fill = "Mutational\nsignature",
color = "Mutational\nsignature") +
ylim(c(-0.5, 1))
plot_box_grid_supplement <- plot_grid(pw_boxplot_patient_jakstat + remove_xaxis,
pw_boxplot_patient_module,
ncol = 1, align = "v",
rel_heights = c(c(0.3, 0.7)))
plot_box_grid_supplement
ggsave("_fig/003_cancer_cell/003_pathway_heatmap_boxplot_patient_lvl.pdf",
plot_box_grid_supplement, width = 8, height = 8)3 session info
devtools::session_info()## ─ Session info ───────────────────────────────────────────────────────────────
## setting value
## version R version 3.6.2 (2019-12-12)
## os Debian GNU/Linux 10 (buster)
## system x86_64, linux-gnu
## ui X11
## language (EN)
## collate en_US.UTF-8
## ctype en_US.UTF-8
## tz Etc/UTC
## date 2021-03-09
##
## ─ Packages ───────────────────────────────────────────────────────────────────
## ! package * version date lib
## abind 1.4-5 2016-07-21 [2]
## ape 5.3 2019-03-17 [2]
## assertthat 0.2.1 2019-03-21 [2]
## backports 1.1.10 2020-09-15 [1]
## bibtex 0.4.2.2 2020-01-02 [2]
## Biobase 2.46.0 2019-10-29 [2]
## BiocGenerics 0.32.0 2019-10-29 [2]
## bitops 1.0-6 2013-08-17 [2]
## broom 0.7.2 2020-10-20 [1]
## callr 3.4.2 2020-02-12 [1]
## car 3.0-8 2020-05-21 [1]
## carData 3.0-4 2020-05-22 [1]
## caTools 1.17.1.4 2020-01-13 [2]
## cellranger 1.1.0 2016-07-27 [2]
## circlize 0.4.10 2020-06-15 [1]
## cli 2.0.2 2020-02-28 [1]
## clue 0.3-57 2019-02-25 [1]
## cluster 2.1.0 2019-06-19 [3]
## codetools 0.2-16 2018-12-24 [3]
## colorblindr * 0.1.0 2020-01-13 [2]
## colorspace * 1.4-2 2019-12-29 [2]
## ComplexHeatmap * 2.2.0 2019-10-29 [1]
## cowplot * 1.0.0 2019-07-11 [2]
## crayon 1.3.4 2017-09-16 [1]
## curl 4.3 2019-12-02 [2]
## data.table 1.12.8 2019-12-09 [2]
## DBI 1.1.0 2019-12-15 [2]
## dbplyr 2.0.0 2020-11-03 [1]
## desc 1.2.0 2018-05-01 [2]
## devtools 2.2.1 2019-09-24 [2]
## digest 0.6.25 2020-02-23 [1]
## dplyr * 1.0.2 2020-08-18 [1]
## ellipsis 0.3.1 2020-05-15 [1]
## evaluate 0.14 2019-05-28 [2]
## fansi 0.4.1 2020-01-08 [2]
## farver 2.0.3 2020-01-16 [1]
## fitdistrplus 1.0-14 2019-01-23 [2]
## forcats * 0.5.0 2020-03-01 [1]
## foreign 0.8-74 2019-12-26 [3]
## fs 1.5.0 2020-07-31 [1]
## future 1.15.1 2019-11-25 [2]
## future.apply 1.4.0 2020-01-07 [2]
## gbRd 0.4-11 2012-10-01 [2]
## gdata 2.18.0 2017-06-06 [2]
## generics 0.0.2 2018-11-29 [2]
## GetoptLong 1.0.2 2020-07-06 [1]
## ggplot2 * 3.3.2 2020-06-19 [1]
## ggpubr * 0.4.0 2020-06-27 [1]
## ggrepel 0.8.1 2019-05-07 [2]
## ggridges 0.5.2 2020-01-12 [2]
## ggsignif 0.6.0 2019-08-08 [1]
## GlobalOptions 0.1.2 2020-06-10 [1]
## globals 0.12.5 2019-12-07 [2]
## glue 1.3.2 2020-03-12 [1]
## gplots 3.0.1.2 2020-01-11 [2]
## gridExtra 2.3 2017-09-09 [2]
## gtable 0.3.0 2019-03-25 [2]
## gtools 3.8.1 2018-06-26 [2]
## haven 2.3.1 2020-06-01 [1]
## hms 0.5.3 2020-01-08 [1]
## htmltools 0.5.1.1 2021-01-22 [1]
## htmlwidgets 1.5.1 2019-10-08 [2]
## httr 1.4.2 2020-07-20 [1]
## ica 1.0-2 2018-05-24 [2]
## igraph 1.2.5 2020-03-19 [1]
## irlba 2.3.3 2019-02-05 [2]
## jsonlite 1.7.1 2020-09-07 [1]
## KernSmooth 2.23-16 2019-10-15 [3]
## knitr 1.26 2019-11-12 [2]
## labeling 0.3 2014-08-23 [2]
## lattice 0.20-38 2018-11-04 [3]
## lazyeval 0.2.2 2019-03-15 [2]
## leiden 0.3.1 2019-07-23 [2]
## lifecycle 0.2.0 2020-03-06 [1]
## listenv 0.8.0 2019-12-05 [2]
## lmtest 0.9-37 2019-04-30 [2]
## lsei 1.2-0 2017-10-23 [2]
## lubridate 1.7.9.2 2020-11-13 [1]
## P magick * 2.4.0 2020-06-23 [?]
## magrittr * 2.0.1 2020-11-17 [1]
## MASS 7.3-51.5 2019-12-20 [3]
## Matrix 1.2-18 2019-11-27 [3]
## memoise 1.1.0 2017-04-21 [2]
## metap 1.2 2019-12-08 [2]
## mnormt 1.5-5 2016-10-15 [2]
## modelr 0.1.8 2020-05-19 [1]
## multcomp 1.4-12 2020-01-10 [2]
## multtest 2.42.0 2019-10-29 [2]
## munsell 0.5.0 2018-06-12 [2]
## mutoss 0.1-12 2017-12-04 [2]
## mvtnorm 1.0-12 2020-01-09 [2]
## nlme 3.1-143 2019-12-10 [3]
## npsurv 0.4-0 2017-10-14 [2]
## numDeriv 2016.8-1.1 2019-06-06 [2]
## openxlsx 4.1.5 2020-05-06 [1]
## pbapply 1.4-2 2019-08-31 [2]
## pillar 1.4.6 2020-07-10 [1]
## pkgbuild 1.0.6 2019-10-09 [2]
## pkgconfig 2.0.3 2019-09-22 [1]
## pkgload 1.0.2 2018-10-29 [2]
## plotly 4.9.1 2019-11-07 [2]
## plotrix 3.7-7 2019-12-05 [2]
## plyr 1.8.5 2019-12-10 [2]
## png 0.1-7 2013-12-03 [2]
## prettyunits 1.1.1 2020-01-24 [1]
## processx 3.4.2 2020-02-09 [1]
## ps 1.3.2 2020-02-13 [1]
## purrr * 0.3.4 2020-04-17 [1]
## R.methodsS3 1.7.1 2016-02-16 [2]
## R.oo 1.23.0 2019-11-03 [2]
## R.utils 2.9.2 2019-12-08 [2]
## R6 2.4.1 2019-11-12 [1]
## RANN 2.6.1 2019-01-08 [2]
## rappdirs 0.3.1 2016-03-28 [2]
## RColorBrewer 1.1-2 2014-12-07 [2]
## Rcpp 1.0.4 2020-03-17 [1]
## RcppAnnoy 0.0.16 2020-03-08 [1]
## RcppParallel 4.4.4 2019-09-27 [2]
## Rdpack 0.11-1 2019-12-14 [2]
## readr * 1.4.0 2020-10-05 [1]
## readxl * 1.3.1 2019-03-13 [2]
## remotes 2.1.0 2019-06-24 [2]
## reprex 0.3.0 2019-05-16 [2]
## reshape2 1.4.3 2017-12-11 [2]
## reticulate 1.14 2019-12-17 [2]
## rio 0.5.16 2018-11-26 [1]
## rjson 0.2.20 2018-06-08 [1]
## rlang 0.4.8 2020-10-08 [1]
## rmarkdown 2.0 2019-12-12 [2]
## ROCR 1.0-7 2015-03-26 [2]
## rprojroot 1.3-2 2018-01-03 [2]
## rstatix 0.6.0 2020-06-18 [1]
## rstudioapi 0.11 2020-02-07 [1]
## rsvd 1.0.3 2020-02-17 [1]
## Rtsne 0.15 2018-11-10 [2]
## rvest 0.3.6 2020-07-25 [1]
## sandwich 2.5-1 2019-04-06 [2]
## scales 1.1.0 2019-11-18 [2]
## sctransform 0.2.1 2019-12-17 [2]
## SDMTools 1.1-221.2 2019-11-30 [2]
## sessioninfo 1.1.1 2018-11-05 [2]
## Seurat * 3.1.2 2019-12-12 [2]
## shape 1.4.4 2018-02-07 [1]
## sn 1.5-4 2019-05-14 [2]
## stringi 1.5.3 2020-09-09 [1]
## stringr * 1.4.0 2019-02-10 [1]
## survival 3.1-8 2019-12-03 [3]
## testthat 2.3.2 2020-03-02 [1]
## TFisher 0.2.0 2018-03-21 [2]
## TH.data 1.0-10 2019-01-21 [2]
## tibble * 3.0.4 2020-10-12 [1]
## tidyr * 1.1.2 2020-08-27 [1]
## tidyselect 1.1.0 2020-05-11 [1]
## tidyverse * 1.3.0 2019-11-21 [2]
## tsne 0.1-3 2016-07-15 [2]
## usethis 1.5.1 2019-07-04 [2]
## uwot 0.1.5 2019-12-04 [2]
## vctrs 0.3.5 2020-11-17 [1]
## viridis * 0.5.1 2018-03-29 [2]
## viridisLite * 0.3.0 2018-02-01 [2]
## withr 2.3.0 2020-09-22 [1]
## xfun 0.12 2020-01-13 [2]
## xml2 1.3.2 2020-04-23 [1]
## yaml 2.2.1 2020-02-01 [1]
## zip 2.0.4 2019-09-01 [1]
## zoo 1.8-7 2020-01-10 [2]
## source
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.3)
## CRAN (R 3.6.2)
## Bioconductor
## Bioconductor
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## Github (clauswilke/colorblindr@1ac3d4d)
## R-Forge (R 3.6.2)
## Bioconductor
## CRAN (R 3.6.2)
## CRAN (R 3.6.3)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.3)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.3)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.3)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.3)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## Bioconductor
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.3)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.3)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.3)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.3)
## CRAN (R 3.6.3)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.3)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
## CRAN (R 3.6.2)
##
## [1] /home/uhlitzf/R/lib
## [2] /usr/local/lib/R/site-library
## [3] /usr/local/lib/R/library
##
## P ── Loaded and on-disk path mismatch.